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Abstract 

Exact analyses are given for two three-dimensional lattice systems: A system 
of close-packed dimers placed in layers of honeycomb lattices and a layered 
triangular-lattice interacting domain wall model, both with nontrivial inter- 
layer interactions,. We show that both models are equivalent to a 5- vertex 
model on the square lattice with interlayer vertex- vertex interactions. Using 
the method of Bethe ansatz, a closed- form expression for the free energy is ob- 
tained and analyzed. We deduce the exact phase diagram and determine the 
nature of the phase transitions as a function of the strength of the interlayer 
interaction. 
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I. INTRODUCTION 



An important milestone in the field of exact solutions of lattice-statistical systems is 
the solution of close-packed dimers on planar lattices obtained by Kasteleyn [jlj and by 
Fisher 0. However, there has since been very little progress in extending the dimer solution 
to higher dimensions. To be sure, Bhattacharjee et al. [[J have studied dimers on a certain 
three-dimensional (3-D) lattice using numerical means, and two of us M have solved a vertex 
model in arbitrary d dimension, a solution which also solves a dimer problem in d dimension. 
In the latter case, however, the dimer model involves unphysical negative statistical weights. 

In a recent Letter ||, hereafter referred to as I, three of us reported on the solution of 
a 3D dimer system as an instance of a more general class of soluble 3D lattice-statistical 
problem. In contradistinction with other exactly solved 3D systems HJT] which invariably 
involve negative Boltzmann weights, the formulation reported in I, which generalizes other 
special cases reported elsewhere ||, marks the success of solving a 3D lattice-statistical 
model with strictly positive Boltzmann weights. In this paper we present details of the 
analysis. In addition, we show also that our solution solves a layered domain wall model 
with interlayer interactions. 

This paper is organized as follows. In Sec. II we define a layered dimer system with 
interlayer interactions and its equivalent layered 5-vertex model. The description of an 
equivalent layered domain wall model is given in Sec. III. The free energy of the 3D system 
is analyzed in Sec. IV with the phase diagram obtained in Sec. V. The critical behavior is 
deduced in Sec. VI. Finally, in Sec. VII, we discuss the occurrence of infinite degeneracy of 
orders in the system. 

II. A LAYERED DIMER SYSTEM AND THE EQUIVALENT 5-VERTEX MODEL 

Consider a 3D lattice £ consisting of K layers of honeycomb lattices stacked together as 
shown in Fig. 1. Each layer of £ is an honeycomb dimer lattice in which dimers with weights 
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u,v,w are placed in the three respective lattice directions. The dimers are close-packed 
within each layer and interact with an interlayer interaction shown in Table I which gives, 
for example, the interaction energy 2h/3, and hence a Boltzmann factor e~ 2h ^ 3 , between a u 
dimer in the kth layer and a v dimer in the (k + l)th layer. This completes the description 
of our 3-D dimer system. Since a perusal of Table I shows that the negation of h corresponds 
to the interchange of the layers k and k + 1, we can without loss of generality take h > 0. 

The honeycomb dimer system can be formulated as a 5-vertex model on a square lattice 
||. This can be seen by drawing the honeycomb lattice in the form of a "brick- wall" as 
shown in Fig. 2. The shrinking of each box containing two lattice points connected by a w 
edge into a point then converts the honeycomb lattice into a square lattice. By regarding 
the presence of a u or v dimer on the remaining honeycomb edges as being bonds, each 
dimer configuration is then mapped into a vertex configuration of a 5-vertex model, and 
vice versa. The resulting five vertex configurations and weights || are shown in Fig. 3. 

Furthermore, the interlayer dimer interaction leads to an interlayer vertex interaction. 
It turns out that the interlayer vertex-vertex interaction corresponding to Table 1 is not 
unique. To deduce a useful interlayer vertex-vertex interaction we first modify Table 1 by 
replacing the uu and vv entries by 2eh, where e = +1 (—1) for sites in sublattice A (B). Since 
two interacting uu or vv dimers are always parallel covering a pair of A and B sites, this 
replacement does not alter the overall interaction energy. A little algebra then shows that the 
dimer interaction of the modified Table I leads to the interlayer vertex interactions shown 
in Table II. Thus, we have at hand a layered 5-vertex model with a particular interlayer 
vertex-vertex interaction. 

Let each square lattice be of size M x N, with M sites in a column and N sites in a 
row. This corresponds to MNK dimers on C. Label sites of the layers of square lattices by 
indices {m,j, k}, with m — 1 . . . M, j = 1 . . . N and k — 1 . . . K. Denote the vertex weight 
at site {m,j,k} by W m jk, and denote the interaction Boltzmann weight between vertices 
{m,j, k} and {m,j, k + 1} as given in Table II by B m ^. Then, it is our goal to evaluate the 
partition function 
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K M N , v 

Zmnk = X! II II II ( BmjkW m jk J (1) 

config. fc=l m=l j=l 

where the summation is taken over all dimer, or vertex, configurations, and the per-dimer 
free energy 

/ = K~ l lim (MN)- 1 lnZ MNK . (2) 

Ai,iV-+oo 

For simplicity, we shall assume K = 3x integers. We also assume periodic boundary condi- 
tions. 

To write the interlayer vertex interactions of Table II in the form of B m jk, we introduce 
variables a m jk = ±1 and (3 m jk = ±1, respectively, for the horizontal and vertical edges 
within the fcth layer and originating from the site {m,j, k} in the direction of, say, decreasing 
{m,j}, such that a m jk = +1 (—1) corresponds to the edge having a bond (empty). It is 
then straightforward to verify that the vertex-vertex interactions in Table II can be written 
as 

£ = -hiptjPj - a j+1 (3j) - ^(aj - a j+1 ) - - (3) 
where we have, for convenience, suppressed the subscripts m and k by adopting the notation 

Pm+l,j,k Pji (3m,j,k+l ~^ Pji (4) 

and similarly for the a's. Now the second and third terms in (|3[) are cancelled upon in- 
troducing this interaction into the overall partition function (|1]). This leads to an effective 
Boltzmann factor 

B mjk = exp(/i(a j /3 j - a j+1 f3j)), (5) 

which is to be used in ([J). 

III. A LAYERED DOMAIN WALL MODEL 

In this section we show that the layered dimer and 5-vertex models of the preceeding 
sections also describe a layered domain wall model with interlayer interactions. 
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Consider a 3D lattice consisting of K layers of triangular lattices whose faces are elemen- 
tary (up- and down-pointing) triangles. Sites of the triangular lattices are occupied by Ising 
spins a = ± with the constraint that, around each face of the lattice, there are precisely two 
spins of the same sign and one spin of the opposite sign. The allowed spin configurations 
are those of the ground state of an isotropic antiferromagnetic Ising model. Furthermore, 
if one erases lattice edges connecting two spins of the same sign, one arrives at a diamond 
(or rhombus) covering of the triangular lattice. This can be interpreted as a dimer covering 
of the dual honeycomb lattice, by placing dimers connecting the two dual lattice points on 
the elongated diagonal of each rhombus. It is clear that the mapping between the spin 
configurations and the diamond and dimer coverings is two to one. Indeed, this mapping 
has been used to extract the solution of the honeycomb dimer lattice from the Ising ground 



state |0|. 



The spin configurations can also be viewed as representing domain wall configurations 
njnj- This mapping is most conveniently seen [II] from the associated diamond covering 



scheme. If one attaches strips to those diamonds oriented in two of the three possible 
directions as shown in Fig. 4, then the strips form continuous lines and propagate in a zigzag 
but generally vertical direction, which can be interpreted as representing domain walls. (Cf. 
Figs. 2 and 4 of for a typical domain wall configuration.) A spin configuration is thus 



mapped into a domain wall configuration. Specifically, the triangular faces of the lattice can 
be in one of the six "strip" configurations shown in Fig. 5, and the domain wall model is 
defined by associating weights to the triangles as shown. 

Next we introduce interlayer domain wall interactions. Shift the (k + l)-th layer by half 
lattice constant to the left with respect to the k-th layer so that the up-pointing (down- 
pointing) triangles in the layer k will be adjacent to down-pointing (up-pointing) triangles 
in the layer k + 1. Let two adjacent triangular faces in planes k and k + 1 interact with 
an energy shown in Table III. Together with the triangle weights given in Fig. 5, this 
completely defines the layered domain wall problem. More precisely, the partition function 
for the domain wall problem is now given by (|l|), where the summation extends over all 



domain wall configurations, with W m jk representing the product of the triangle weights 
given in Fig. 5 and B m jj. the interlayer interaction given by Table III. 

The mapping of domain wall configurations to five-vertex arrow configurations has been 



given in ||11|| , where the triangular lattice was deformed into a square lattice by tilting it 
clockwise, leading to a 5- vertex model with co> 3 = (instead of u\ = as in Fig. 3). For the 
present paper, we deformed the triangular lattice by tilting it counterclockwise. Then, the 
vertex weights reduce exactly to those given in Fig. 3. 

To obtain an explicit form for B m jk, it is straightforward to verify that, in the language 
of layered 5- vertex model, the interaction of Table III can be written as 

e = -hfajftj - aj+iPj) + h{aj - a j+1 ) - fo(/3j - /%). (6) 

Again the second and third terms in ([]) are cancelled in the overall partition function ([[]), 
so that the effective interaction Boltzmann factor now assumes the form 

B mjk = exp^ajP'j - a j+1 /3j)), (7) 

which differs slightly from ([5]) for the dimer problem. However, repeating precisely the 
same line of argument as in I, one can show that the interlayer interaction (0) leads to 
precisely the same free energy (||) and (||) given below. Thus, the domain wall problem 
(with interlayer interactions of Table III) is completely equivalent to the dimer system (with 
interlayer interactions of Table I). 



IV. THE FREE ENERGY 

In the preceeding sections we have established the complete equivalence of the layered 
dimer and domain wall problems, and their further equivalence with a layered 5-vertex 
model. In this section we analyze the free energy of the layered 5-vertex problem. For 
simplicity we use the language of the dimer system. 

It has been shown in I that the layers of 5-vertex models with interlayer interaction ([S]) 
can be solved by applying a transfer matrix in the vertical direction and a global Bethe 
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ansatz consisting of the usual Bethe ansatz within each layer. This leads to the following 
expression for the free energy 

f(u,v,w,h)= max f({y k }), (8) 
-i<i/fc<i 

where 

1 K 1 /■7r(l-l/ fe )/2 
K k=1 TT JO 

Here, 

■y N y N 

3=1 3=1 

is a quantity conserved from row to row (of vertical edges) in the kth layer square lattice. 
Specifically, we have y^ = 1 — 2n^/N, where is the number of vacant edges in a row. 
Analysis leading to (|9[) has been given in I and will not be reproduced here. 

It is clear that for large u, v or w, the system is frozen with complete ordering of u, v, 
or w dimers in all layers, and hence the free energies 

fu = hiu, U phase 

fv = hit?, V phase 

fw = lnw, W phase. (10) 

These are frozen orderings which we refer to as the U, V, and W phases, respectively. For 
large h, it is readily seen from Table 1 that the energetically preferred state is the one in 
which each layer is occupied by one kind of dimers, u, v, or w, and that the layers are 
ordered in the sequence of {u, w, v, u, w, v ■ ■ •}. This ordered phase is referred to as the H 
phase with the free energy 

f H = -ln(uvwe* h ), H phase (11) 

3 

obtained from a perusal of Table 1. 

For any layer with = +1 the corresponding integral in (Q) vanishes, and for y^ = — 1 
the integral can be evaluated using the integration formula 

7 



w v 
- + -e- 

u u 



d9. 



(9) 



- r ln\ A + Be i9 \d6 = max{m \A\, In 

TT JO 



(12) 



Therefore, the free energy (^) can be explicitly evaluated if y& = ±1 for all fc. Further 
discussion of this case will be given in Sec. VII. 

We have carried out analytic as well numerical analyses of the free energy (§) for fixed 
u,v,w, h, and it was found that the set {yk} which gives the extremum value in (H) always 
repeats in multiples of 3, namely, satisfying fl3]| 



Vk+3 = Vk- 



The following extremum sets of {yk} are found: 

1- 0/i, 2/2, 2/3} = {1,1, 1} : In this case we have all yk = 1, and hence from 



/ = /« 



U phase 



(13) 



This gives rise to the U phase. 

2. {1/1,2/2,2/3} = { — 1, — 1, —1}: In this case we have all = —1, and hence from 



/ = In u H — / In 

7T JO 

= fw, 
= fv, 



W V 

- + -e 

u u 



W > V 

V > w. 



de 



W phase 
V phase 



(14) 



This gives rise to the W and V phases. 

3. {2/1, 2/2j 2/3} — {1, — 1, — 1}: Substituting this sequence of y^ values into (^) and making 
use of (P^D in the resulting expression, one obtains 



/ = lnu + — / In 

D7T J -tv 
= T^fu + 3/^' 

= f H . 



w v 



_ e 2h e ie 



u u 



d9 + — I In 

07T 



-2hit 



W V 

— + -e _i "e" 
u u 



de 



ve 4h < ve 4h < w 



—\h Ah 

w <ve < ve 



ve 4h < w < ve 4h , H phase 



(15) 
(16) 
(17) 



Now the free energies flT5|) and (|T^) can be discarded since they are always smaller than the 
largest of {fu, fv, fw}- Thus, this set of {yk} leads to a frozen ordering for sufficiently large 
h as indicated in (|TTD, which is the H phase. 



4. {2/1, 2/2? 2/3} = {y,y,y}'- In this case all yk = y, where y maximizes the free energy 
(§). Then, substituting yk = y into (H) and carrying out the maximization in (||) by a 
straightforward differentiation with respect to y, one obtains 





— + -e lf> 


/ In 


/o 


u u 



Y phase 



(18) 



where the extremum y is given by 
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[1 - y ) = cos 



-1 



222 
11 —w —v 



2wv 



(19) 



This is a disorder phase which we refer to as the Y phase. Despite its apparent asymmetric 
appearance, the free energy fy(y) is actually symmetric in u,v,w. Note that, for large 



v ~ u>, we have 



71 

^-(1 — 2/o) = k-9 , 



(20) 



where 8 Q is small and given by 



(21) 



5. {yi, 2/2, 2/3} = {2/1, — 1, — 1} : This is the phase with the u layers replaced by layers with 
Vk — so that the layer ordering is {yi, w, v, yi, w,v ■ ■ •}. This is a partially ordered phase 
which we refer to as the I u phase. Again, the substitution of these values of {yk} into (0) 
and a straightforward maximization yield, after using 



/ = fi u {yio) = 3 

_ 1 

" 3 
1 

" 3 



2fw + fy(yio)\, w>ve 4h ^ 1+ y^ 
2fv + fr(yio)], w<ve- 2h ^ 1+ y^ 



ve 



-2/1(1+1/10) 



< w < v e 



,2/1(1+3/10) 



I u phase 



(22) 
(23) 

(24) 



where fy (yio) is defined in (|18"D , with the extremum y w given by 
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i - yw) = cos 1 



u 2 e 8h — w 2 — f 2 " 1 



2wv 



(25) 



The free energies ( p2|) and (^3|) are discarded since they are always less than the largest of 
{fw, fv, fy(yw)}, and we have fy(y w ) < fy(yo) by definition. Therefore, the free energy of 
the I u phase is given by (|2l|). Note that, for large v ~ w, we have 

l(l-y w ) = n-9 1 , (26) 

where Q\ is small and given by 

Q\ = [ u 2 e 8h -(w- v) 2 }/wv. (27) 

5. {yi, y%, 2/3} = {1,2/2, — 1} ; This is the if phase with the w layers replaced by layers with 
2/fc+i — 2/2, so that the layer ordering is {u, y 2 , v, u, y 2 , v ■ ■ ■}. We refer to this as the I v phase. 
Due to the intrinsic symmetry of the interlayer interaction, the free energy of the I w phase 
is the same as that of I u , given in (|23|), with the cyclic permutation of u — > w — > v — > u. 
Alternately, one can substitute these values into (|9]) and carry out the maximization. 
It can be verified that this leads to 

1 



/ = fi w (2/20) 



2fu + fv + 



/■t(1-J/2o)/2 


v + 




/ In 


dO 


Jo 


u u 





I w phase 



with 
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(1 - 2/20) = cos 1 



u 2 — w 2 e 8h — v 2-1 



(28) 



(29) 



2wve ih 

6. {1/1,2/2,2/3} = {1, — 1,2/3}: This is the H phase with the v layers replaced by layers with 
Vk+2 = 2/3, so that the layer ordering is {u, w, y 3 , u, w,y 3 - ■ •}. We refer to this as the I v phase. 
Again, the free energy of the I v phase can be written down by symmetry. Alternately, the 
substitution of these values of {y^} into @ and (H) yields 

1 



/ = //„(2/3o) 



2fu + fw + 



/•7r(l-2/3o)/2 






/ In 




/o 


u u 





I v phase 



with 



7T 



(1 - y 30 ) = cos 1 



u 2 — w 2 — v 2 e 8hl 



2wve 4h 



(30) 



(31) 
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V. THE PHASE DIAGRAM 



Since the phase diagram must reflect the {u, v,w} symmetry of the interlayer interaction 
given in Table 1, it is convenient to introduce coordinates 

X = \n(v/w) Y = (Vzy 1 \n(vw/u 2 ) (32) 

such that any interchange of the three variables u, v, and w corresponds to a 120° rotation in 
the {X, Y} plane. The phase boundaries are then determined by equating the free energies 
of adjacent phases. The results are collected in Fig. 6. 

The phase diagrams depend on the value of h and are different in different regimes. 

• h < h : For small h the phase diagram is the same as that of the h = noninteracting 
2D system, namely, the diagram shown in Fig. 6a. The phase boundary between the 
{U, V, W} phases and the Y phase, which stays the same in all regimes below, is obtained 
by setting y = ±1 in fllTf ) where fu, fy or f w is equal to fy. These boundaries are 

u = \v ±w\. (33) 

In terms of the coordinates X and Y, u = v + w and u — \v — w\ read, respectively, 

Y = ^ln[2cosh(X/2)], Y = -| ln[2 sinh(|X|/2)]. 

• ho < h < hi. As h increases from zero, our numerical analyses indicate that the H 
phase appears when h reaches a certain value h Q . The resulting phase diagram is shown 
in Fig. 6b. The phase boundary between the H and Y phases is given by /# = fy, or, 
explicitly 

l -\n{uvwe ih ) = f Y (y ). (34) 

Thus, ho is obtained from (|34"D by setting u = v = w (X = Y = 0) where the H phase first 
appears. This yields 7r(l — yo)/2 = 27r/3 and 

3 f 2w / 3 

h = — ln(2 + 2 cos 9)d9 = 0.2422995 • ■ ■ . (35) 

07T JO 
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• h± < h < hi'- As h increases from ho, the I U ,I V ,I W phases appears when h reaches a 
certain value hi. The resulting phase diagram is shown in Fig. 6c. Now the I u phase is the 
H phase with the u layers (with y^ = 1) replaced by layers with y^ = y i, the boundary 
between the two regimes is therefore given by yio = 1 or, explicitly using (p5|), 



Ah 

w + v = ue . 

The boundary between the I u and the Y phases is fyiyo) = fi u (yio) or , explicitly 

1 Ml-Vo)/2 



(36) 



In it 



7T JO 



In 

2/1(1+3/10) 



v w 
- + —e 
u u 



,0 



dd 



v 1 /■7r(l-j/lo)/2 


v w w 


de 


+ In u H — / In 


/ n Jo 


u u 





(37) 



The boundaries of the I v and I w regimes can be written down similarly. 

To compute the numerical value of hi, we note that the I u phase first appears at v = w 
(X = 0) when all three phases H, Y and I u coincide. Therefore, hi is obtained by solving 
(HP and (]36[) for v/u and h at w — v. This leads to the value 



h 



= 0.2552479 • • • 



where v/u is the solution of the equation 



l+j/o ln- + -(1 + 3j/o ln2 = - / 
■u 6 n Jo 



ln(l + cos 



(38) 



with 



7T 



;i - 2/0) 



cos 



-1 



u 



1 



(39) 



2 V " " v/ V2w 2 
• h 2 < h < h 3 : As h increases from hi, it was found that the regimes I u , I v , and I w 
extends to infinite when h exceeds a certain value h 2 . The phase diagram is shown in Fig. 
6d. The value of h% can be deduced from (pTj) in its large w = v expansion. Setting w = v 
in (|37D, introducing ( |20| ) and ( p6|) for large w, v, and equating the coefficients of \n(v/u) on 
both sides of the equation, one obtains 



1 



7T 



7T 



2 + -(tt-0i) 

7T 



(40) 



12 



or, simply, 6 Q = 9i/3. Now from (|21|) , fl2~7|) , and w = t> , we have the expressions #0 = u/v 
and 6 1 ! = (u/v)e 4h . It follows that we have 

ha = (In 3)/4 = 0.3816955- ••. 

• h 3 < h < h±: As h increases from h 2 , it was found that the boundary of the H 
phase bulges toward the U, V, W phases along the 30°, 150°, 270° lines, touching the U, V, W 
boundaries in these directions when h equals a certain value h^. For h > h 3 , the H phase 
borders directly with the U, V, W phases with respective boundaries 



u 2 = vwe 4,11 , 



v 2 = wue 4h , 



w 2 = uve 4h . 



(41) 



The size of these borders grows while the Y phase shrinks as h increases. The phase diagram 
in this regime is shown in Fig. 6e. To determine h 3 , we let the {H, Y} phase boundary fl34|) 
to touch the {Y, U} phase boundary u = w + vaXw = v (the 270° direction). Using ( [19]) 
we have y = 1, and it follows that (Rl) becomes 



\n(2v d e 



3 4h\ 



ln(2v) + 



from which we find 14 



h 3 = (ln2)/2 = 0.3465735 ••• . 

• h > h^. As h increases further from h 3 , it was found that the Y phase disappears 
completely when h exceeds a certain value h^. The phase diagram in this regime is shown 
in Fig. 6f. To determine the numerical value of h^, we note that the Y phase disappears 
when the boundary v = w + u between the V and Y phases coincides with the boundary 



( |37|) between the I u and Y phases at large w,v. Therefore, we again expand (|37D for large 
v,w but now subject to v — w = u. Introduce ( p0|) and fl26|) for the integration limits. But 
now from (^l|) and ([27D we have 9q = 0, d\ = ju/ \/wv, where 



1 = ye 8h -l. (42) 
Substituting ( pO"D and ( p^) into d37|) and making use of (|12| ) and the relation (for v > w) 



1 


rir—01 


V 


W if) 




/ In 




+ -e id 


7T 


Jo 


V 


u 



V 


1 


/ In 


V 




II 


7T 


Jo 


u 


u 
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one obtains 



Inv = - 

3 



, / Aha fir\ , 1 f dl , f v 2 w 2 2vw „\ 

ln(vwe 4Ml/ ^) + In w / In — + — cos 9 d6 

2ix Jo \u 2 u 2 u 1 J 



(43) 



Since Q\ is small, we can write cos# = 1 — 9 2 /2 in the integrand, and the integral can be 
simplified by introducing the change of variable y = \/wv9/u. Introduce w = v — u and 
expand (|43] ) for large w, v using, for example, \n(vw) = 2 In v — u/v, the leading terms of the 
order of hav are cancelled. The next terms including the integral are of the order of 0(u/v). 
Setting the coefficient of these term equal to zero, one obtains 



4/i 7 = vr + i r\n(l + y 2 )dy 
2 Jo 



= ^ + ^7 ln(l + 7 2 ) - 7 + tan 1 7, 

or, after using 1 + r y 2 = e 8h , 

7 — tan -1 7 = 7r (44) 
whose solution gives h^. Specifically, we find 

h 4 = - ln(l + 7 2 ) = 0.3816955 • • ■ . 
8 

Phase diagram for the domain wall model. 

Since the domain wall model with weights given in Fig. 5 and interlayer interactions 
of Table III is completely equivalent to the layered dimer system, the phase diagram of 
the domain wall model is the same as that of the dimer model. For example, the phase U 
corresponds to a phase with no domain walls, and that the phase H corresponds to a sequence 
of triplets of layers with no domain walls, maximal density of domain walls consisting of 
elementary weights y/w, and maximal density of domain walls of weights \/v (Cf. Fig. 5). 



VI. THE CRITICAL BEHAVIOR 



In this section we determine the critical behavior near all phase boundaries. Since the free 
energies are given by different analytic expressions in different phase regimes, one generally 
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expects the first derivatives of the free energy (with respect to a temperature T, say) be 
discontinuous. This then leads to first-order transitions. However, if the first derivatives of 
the free energies happen to vanish on both sides of the boundary, then one has a continuous 
transition. Applying this reasoning, we find all transitions to be of first order, except those 
between the {U, V, W} and Y phases, and between the {I u , I v , I w } and H phases, which are 
found to be the same as the transition in the 5- vertex model |9|,[15|, namely, a continuous 
transition with a square- root singularity in the specific heat. This transition, first reported 



by one of us in 1967 [fL5|| , is now known as the Pokrovsky-Talapov type phase transition ||16|| . 

Regarding u, v, w and e h as Boltzmann factors, the ordered U, V, W and H phases (with 
y& = ±1) have constant free energies and hence zero first derivatives. Therefore, we focus 
on the boundaries of these frozen regimes. 

We have seen that the transition between the U, V, W phases and the Y phase is the 
same as that occurring in the 2D system, which is known || to be continuous. This fact 
can also be seen by expanding the free energy near y Q as 

fy(y) = fv(yo) + (y- y )f Y (y ) + ^(y - yo) 2 f Y (y ) + ~(y - yo) 3 f Y ' (y ) + • • ■ ■ (45) 



Using the expression of fy(y) defined by (p~8|) , one sees that, indeed, the first derivative 
f Y (yo), yo — ±1) vanishes identically on the boundary (p3|) which is precisely f Y {yo) = 0. 
Furthermore, it is also seen that fy{yo) ~ sin[7r(l — yo)/2] = 0. Therefore, the extremum of 
f Y (y) (H) occurs at y = y cxtrm given by 



2/extrm — 2/0 — it 



\-fv'(yo) ' 1 j 



where t — \T — T c \, T c being the critical temperature. Substituting this y ex trm into 
obtains 

2,/ 



<Y 



2/y(2/o) 



\ -fy'iyo) 

= fy(y )+c(u,v,w,h)t 3 / 2 . (47) 

This leads to a square-root singularity in the specific heat, which is a characteristic of the 
Pokrovsky-Talapov phase transition. The key element leading to this result is the fact that 
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the boundary of the frozen phases is given precisely by fy{yo) = 0, rendering the first 
derivative of the free energy to vanish at the boundary. 

Applying the same analysis to the H and I u phases, the boundary is again given by 
fi u (yio) = 0. Furthermore, it is also seen that f" u (y 10) = identically It follows that the 
analysis can be carried through exactly as given in the above, and one concludes 

//.G/extan) = /(yio) + C!(U,V,W, h)t 3/2 . (48) 

This gives rise again to a square- root singularity in the specific heat. The consideration of 
the I v , I w and the H boundaries can be carried out in a similar fashion. 

VII. DEGENERACY OF ORDERED STATES 

We discuss in this section the degeneracy of ordered states. Particularly, we show that the 
system has a nonzero per-layer entropy on the boundaries between H and U, V, W phases. 

We first establish the occurrence of a degeneracy from an energy consideration. For this 
purpose it is sufficient to show that this is the case along the {H, V} boundary fl4T|), namely, 

e ih = v 2 /wu. (49) 

We already know that, along the boundary fl49|) , the following layer orderings of the H and 
V phase are degenerate, 

...222222222... 

...(132)(132)(132).... (50) 

Here, for convenience, we have used the notations {1,2,3} for {u,v,w}. Generally, when 
y k = ±1, each layer contains dimers of only one kind, u,v or w. Let a iy i — 1,2,3 denote 
the numbers of u, v and w layers, respectively, as a fraction of K, and a%2 the fraction of 
adjacent pair of u, v layers, etc. Then, perusal of Table I shows that this leads to the per-site 
free energy 
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f = \ f «i In u + a 2 In v + a 3 In w] + ^ ( a 13 + a 32 + a 2 i - a 12 - "23 - "31 ) • (51) 



2V J 3 

Here, the a's satisfy the conservation rules a% + a.i + a 3 = 1, ^ a^- = «j. 
Consider the following ordering, 

...2(1)3222(1)322222(1)32..., (52) 

characterized by single u layers separating strings of layers of the type wvvvvv... where there 
is at least one v layer in each string. It is readily seen using fl52"D that for this ordering we 
have 

"12 = "23 = "31 = 

"21 — "32 = "13 = "l = "3 = "> "2 = 1 — 2tt. (53) 

It is then a simple matter to substitute (pE§D, (p^) into (§TJ), obtaining f = f v . Thus, any 
layer ordering having the structure of (^2j) is degenerate to fv on the {V, H} boundary. 
Obviously, the number of such structures is infinite as K — > oo. Note that ( |5T1D is a special 
case of (0). This degeneracy has also been confirmed in our numerical analysis of @. 

Generally ii yk = ±1 for all layers, the free energy (H) can be explicitly evaluated for 
any ordering. Let p a ^', a — ±, denote the fraction of layers with y k = —1 such that 
{yk-i,yk,Vk+i} = W, — , 0"'}, where for brevity we denote ±1 by ±. Consider, for example, 
the domain v < w. A straightforward evaluation of (||) leads to the expression 

/ = a + In u + a_ In w, t> e~ 4/l < f e 4/l < w 

= a + In w + (p ++ + + p ) lnv 

In io + 4hp_ + ve~ ih <w< ve 4h . (54) 

The first line can be discarded since it is always smaller than the larger of Inu and lnw. 
A degeneracy of states now occurs if the second line coincides with the free energy of any 



phase. In the case of the ordering (p2[) , for example, it is translated to 

•••-(+) (+) (+) ■ 
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Now since each + layer is associated with precisely one H — and one — h neighbors, we have 
a + = p + _ = P-+. The substitution of a + = = p_ + and ( |4"9]) into (Q) now leads to 
/ = fv, in agreement with the energy consideration. The degenerate states on the {H, U} 
and {H, W} boundaries can be obtained by cyclic permutations of u, v, w. 

The entropy of the ordered state Q52] ) can be computed. We note that the main feature 
of ( [52]) is that layers v are followed by either u or v layers, and u layers are followed by only 
w layers, and w layers by v. Then the degeneracy S of the sequence (|52|) is given by the 
trace of a transfer matrix as 



S = Tr (T K ) = Af + Xf + Xf ~ \( 



.A' 



K i \K 



K large 



(55) 



where 



T 



/0 1 \ 

1 1 
VO 1 OJ 



and A,'s are the eigenvalues of T with Xi > \XA, j = 2, 3. We find 



,29 3V93\ l/3 {29 3 V / 93\ 1/3 

1 + , t + — ) + (y"^ 



(56) 



1.46557... 



VIII. SUMMARY 

We have considered a three-dimensional layered dimer system with interlayer interac- 
tions and its equivalent layered domain wall model, and analyzed its exact solution. It is 
found that the phase diagram, shown in Fig. 6, depends crucially on the strength of the 
interlayer interactions. There exist ordered U, V, W, H phases corresponding to, respectively, 
large dimer weights u, v, w and large interlayer interactions h. In addition, there also ex- 
ist a disorder phase Y and partially ordered phases I U ,I V ,I W . The phase boundaries are 
determined by equating the free energies of adjacent regimes. Particularly, the boundary 
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between the U, V, W phases and the Y phase assume the simple form (|33|) , between the H 
and Y phases the form fl34]), and between the H and I u phases the form (p6|). All transitions 
are found to be of first order, except the transitions between the U, V, W phases and the Y 
phase, and the transitions between the I U ,I V ,I W and H phases, which are of second order 
with a square-root singularity in the specific heat. 
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FIGURES 



FIG. 1. A three-dimensional lattice model consisting of layers of honeycomb dimer lattices. 
FIG. 2. The mapping of a honeycomb lattice onto a square lattice. 

FIG. 3. Vertex configurations and weights of the 5-vertex model. 

FIG. 4. The three possible orientations of a diamond. Strips are associated with diamonds 
oriented in two particular directions. 

FIG. 5. The six strip configurations and the associated weights of a triangle. 

FIG. 6. Phase diagram of the 3D system. 
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TABLES 

TABLE I. Interaction energy between two dimers incident at the same site of adjacent layers. 
The interaction is symmetric in u, v, w. 

layer k — > k + 1 u v w 

u 2/i/3 -2/i/3 

v -2/i/3 2ty3 

w 2/i/3 -2/i/3 



TABLE II. Interaction energy between two vertex configurations of adjacent layers. 
tOi, i = 2, • • • , 6 denotes the vertex of type i in Fig. 3. 

layer k — k + 1 u;2 ^3 ^4 ^5 ^6 

w 2 -4/i/3 4/i/3 

u 3 Ah/3 -4/t/3 4/t/3 -8/i/3 

w 4 -4/i/3 4/i/3 -4/i/3 8/i/3 

L,j 5 8/i/3 -8/i/3 

w 6 -4/i/3 4/i/3 



TABLE III. Interaction energy between two strip triangles of adjacent layers. The triangle 

configurations are as numbered in Fig. 5. 

layer k -► k + 1 1 2 3 4 5 6 

1 -2h 

2 2h 

3 

4 -2h 

5 2h 

6 
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